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Abstract 

A method is presented for solving elastodynamic problems in radially inhomogeneous elas- 
tic materials with spherical anisotropy, i.e. materials such that Cijki = Cijkiir) in a spherical 
coordinate system {r,6,(p}. The time harmonic displacement field u(r, 0, 0) is expanded in a 
separation of variables form with dependence on 9, cj) described by vector spherical harmonics 
with r-dependcnt amplitudes. It is proved that such separation of variables solution is generally 
possible only if the spherical anisotropy is restricted to transverse isotropy with the principal 
axis in the radial direction, in which case the amplitudes are determined by a first-order ordi- 
nary differential system. Restricted forms of the displacement field, such as u{r,9), admit this 
type of separation of variables solutions for certain lower material symmetries. These results 
extend the Stroh formalism of elastodynamics in rectangular and cylindrical systems to spherical 
coordinates. 

1 Introduction 

The Stroh formahsm [1], which recasts equations of time-harmonic elastodynamics in the form of 
a first-order ordinary diff'erential system (ODS), is a powerful technique for dealing with elastic 
materials inhomogeneous in one coordinate. The method was originally established for rectangular 
coordinates, e.g. [21 [3l [H [5], and has been developed for applications in cylindrical coordinate 
systems [6l[7]. One complicating factor for cylindrical, as compared with rectangular, anisotropy is 
that the radial and azimuthal basis vectors and eg depend upon the angular coordinate 6. This 
however does not hamper separation of variables and allows for a Stroh-like ODS provided that 
the material coefficients depend either on the radial or axial coordinate r or z. The situation is 
quite different for spherical anisotropy. The lowest anisotropy that supports general displacement 
fields which can be described by a separation of variables appears to be transverse isotropy (TI) 
with the axis of symmetry in the radial direction, as shown by [8] for statics and [9] for dynamics 
((the term "spherically isotropic" used by [8] and others is equivalent to TI). A state space system 
was developed for radially inhomogeneous TI by [9] , who also identified two distinct types of wave 
motion solutions: an uncoupled pure shear motion and a coupled radial- angular solution pair. The 
state vector approach has been applied to vibrations of thick-walled TI shells [101 [TT] and further 
developed for piezoelectric shells [12]. Lower symmetry can support specific types of kinematically 



1 



restricted deformation. The most general form of spherical anisotropy which admits static solutions 
of the form u = u{r)er in spherical coordinates (r, 6, <j)) is described in jl3j . 

The purpose of this paper is to present a method for solving elastodynamic problems in radially 
inhomogeneous elastic materials with spherical anisotropy, i.e. materials such that Cijki = Cijki{r) 
in a spherical coordinate system [14] . No a priori restrictions are made on the form of the dis- 
placement field. The main departure from previous works is the use of vector spherical harmonic 
functions as the set of vector basis functions. We show that the most general type of spherical 
symmetry, for which the basis of vector spherical harmonics always yields a separable solution, is 
transverse isotropy about (certain lower symmetries admit such solutions as well but only un- 
der appropriate kinematic restrictions). Such anisotropy restriction may actually be not so severe, 
since any spherical anisotropic material with properties independent of the polar-axis orientation 
must possess transverse isotropy anyway (see below). The spherical TI problem reduces to an ODS 
in the radial coordinate r with a system matrix possessing Hermiticity properties that guarantee 
physical attributes such as energy conservation. The key feature of the present analysis is the set 
of basis functions, vector spherical harmonics, which allow for the first time application of the full 
Stroh formalism to spherical elasticity for arbitrary displacement fields. 

The paper is organized as follows. The concept of spherical anisotropy is revisited in ^ Vector 
spherical harmonic basis functions are introduced in ^ The main result, which is the Stroh- 
like ODS for the state vector comprising the r-dependent components of displacement and radial 
traction in the basis of vector spherical harmonics, is described in 21 ^-iid proved in detail in ^ 
Explicit solutions and their properties are discussed in ^ Conclusions and further prospects are 
presented in ^ 

2 Elastic anisotropy in cylindrical and spherical coordinates 

The concept of cylindrical and spherical elastic anisotropy was introduced by Saint- Venant and 
subsequently developed by Lekhnitskii [14j . It is motivated by the existence of special materials 
possessing either a physically distinguished direction aligned with the axis Z of the cylindrical co- 
ordinate system or a point which can be identified as the origin O of the spherical system. It is 
physically relevant to consider such materials in terms of tensor fields in an orthogonal curvilinear 
coordinate system. A tensor $ of order p G N is defined at every point r of a material body V by the 
array {ij = 1, 2, 3; j = 1, ...,p) of its components in some frame of orthogonal basis vectors 

(ei, e2, 63) associated with r G M^. The components change in the usual manner under a change of 
basis (see e.g. |151 §2.6]). Recall that the frame {e} = (e^,, e^, e^) of a rectangular (orthogonal recti- 
linear) coordinate system in is independent of r and of the origin point O, whereas an orthogonal 
curvilinear coordinate system implies a varying frame, e.g., the frame {e(0)} = (6^,651,62) with a 
fixed longitudinal axis Z \\ for the cylindrical system, and the frame {e (0, 0)} = (e,., eg, e,^) with 
a fixed origin O and polar axis Z for the spherical system. The essential difference between the 
rectilinear and curvilinear arrays of components is related to the meaning of a uniform tensor, i.e., 
one independent of r. If the array of components {*&ii..ip} in some rectangular frame is uniform: 
^i-i..ip =const. for all ij,j and r G V, then the tensor $ is also uniform, i.e. is the same for all r, 
and vice versa. By contrast, a uniform array of components {^ii.Ap} referred to a curvilinear frame 
{e} describes a tensor $ which generally depends upon r since {e} does. For example, $ = $ (0) 
for a uniform cylindrical array of components and $ = $ (0, (/>) for uniform spherical components, 
unless these tensors are isotropic (components invariant to change of basis under SO (3)) or if the 
cylindrical one is transversely isotropic (invariant under 5*0(2) about e^). As a matter of defini- 
tion, a tensor $ associated with cylindrical or spherical components which is not isotropic in the 
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above sense, is called cylindrically or spherically anisotropic, respectively. 

The components of a tensor may certainly be non-uniform so that {^ii.Ap} = {^ii.Ap (r)} (where 
r does not have to be related to the same coordinate system as the frame {e}). Many applications 
deal with a specific type of inhomogeneous body V in which the curvilinear components $ji..ip(r) 
of a non-uniform tensor vary with position but maintain the same anisotropy (symmetry class) at 
every r G V. This is the general framework considered in this paper. A widely encountered example 
is the case of a cylindrical or spherical radially inhomogeneous elastic tensor with Cij^i = cijki (r), 
where r may be measured from any point of the fixed axis Z of the cylindrical system (since c^ji^i 
do not depend on z) or from the fixed origin O of the spherical system. 

Another particular aspect of cylindrical and spherical tensor components {$ii..ip} stems from 
the fact that the orientation of frame vectors is undefined at certain points, namely, (i) of the 
pair (6^,651) of the cylindrical frame at the points r = lying on the longitudinal axis Z, (ii) of 
the pair [eQ,e^) of the spherical frame at the points ^ = 0, vr of the polar axis Z, and (iii) of the 
whole spherical frame (6^,60,6,^) at the point r = of the origin O. Now suppose that a given 
body V contains either the axis Z or the point O and has a tensor $ described by the array of 
components {$j^,,ip(r)} which is single- valued at every r £ V including the above-mentioned special 
points. Then in case (i), the array that is independent of r must at every r be invariant 

to the orientation of the pair (e^je^) , i.e. $ must be transversely isotropic about Z \\ e^; in case 
(ii), {^ii.Ap} that is independent of 6 must at every r be invariant to the orientation of the pair 
{e0,e^), i.e. $ must be transversely isotropic about e^; in case (iii), {$ii..ip} that is independent of 
r must at every r be invariant to the orientation of the spherical frame (e^., e^i, e^), i.e. $ must be 
isotropic. These restrictions on cylindrical or spherical anisotropy, which result from eliminating a 
singularity that exists only at isolated (axial or origin) points, are neither immanent nor physically 
binding and can certainly be circumvented formally, say, by assuming a small cavity surrounding 
each point. However, in the spherical case it is physically reasonable to single out the class of 
spherically anisotropic materials invariant with respect to any orientation of the polar axis which 
may thus be called materials with complete spherical anisotropy. The remedy for case (ii) described 
above must then be provided at every r. Hence, material tensors in a body with complete spherical 
anisotropy can only be either uniform or radially inhomogeneous and, unless isotropic, they must 
be transversely isotropic about e^. They are well defined everywhere except a single origin point 
O where (iii) needs to be addressed. 

3 Governing equations 
3.1 Elastodynamic equations 

The dynamic equilibrium vector equation for a linearly elastic material when expressed in spherical 
coordinates is 

r~'^{rhr),r + {rsmey'^[{sm6t0)^e+t^,^ + sm6Ktg + m^] = pu (1) 
/O -1 0\ 

with K= 1 {=-K'), H= I -cos^l f=-H' 
\0 0/ 

Here p = /9(x) is the mass density, u = u(x, t) the displacement and tj = tj(x, t) = eja {i = r,9,(j)) 
the traction vectors defined by the elements of stress (T{x,t) in the orthonormal basis (6^,651, e^) of 
the spherical coordinates (r, 6, (j)). The left member in ([T]), divo" = Yli V- (ej®tj), follows using Qi 
plus the identities || = a^g + Ka, || = a_<^ + Ha, where the comma suffix notation indicates partial 
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differentiation of components only: a^^ = f = 9,<j). The stress elements in the basis 

(er,e£),e^) are aij = CijkiSki where Cijki = Cijki{^) are the components of the spherically anisotropic 
elastic tensor, eu are the components of the strain e(x,t) = ^ [Vu + (Vu)"^] , with summation on 
repeated indices assumed and T for transpose. The traction vectors can therefore be written as 
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where, in the notation of [TB], the matrix (a6) has components (a&)jfc = aiCijuih for arbitrary 
vectors a and b. 



3.2 Vector spherical harmonics 



Our objective is to develop separation of variables vector solutions in the form v(r, 0, (p) = VA{r)A{9, i 
where the three vectors A (e^) are independent of r and provide a complete basis for representing 
vectorial functions of the spherical angles. Vector spherical harmonics are one such set of functions. 

It is useful to first introduce the angular parts D and D • D = of the vector differential 
operators V and V • V = (= A) in spherical coordinates: 
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^D^{ for pure angular functions /(e^) and f(er). Another 



ingredient is the set of (scalar) spherical harmonics of polar order n and azimuthal order m, 
for which there are several slightly different notations in use. Following \17\ p. 64], let 
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where P^{cos9) are the associated Legendre functions of the first kind. The functions Y^ 
satisfy the equation 
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In these terms, the vector spherical harmonics are defined as 

Bmn{er) = X"^DY;^{er), 



with A = [n(n + 1)]^/^ (5) 
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(see \17\ §3.17] for a literature review). The vector harmonics are pointwise orthogonal 

Pmn ■ ^mn — ^mn ' ^mn — ^mn ' -Pmn — 0, (7) 

and orthonormal when integrated by d 17 = sin 9(19d(j): 

d$7Pmn-B*^= / dilBm„-C*,^= / d 17 Cmn ' P^i/ = 0, (8) 

n JO. 

d 17 Pmn ■ Pfj;/ — / d 17 • B^j^ — / d 17 Cmn ' — ^m^^nuj 

Jn Jn 

where the latter two identities hold at n, 7^ 0, see ^6l !|6.4p . Using ||T8| p. 1901], the vector spherical 
harmonics can be shown to satisfy the following identities where / = /(r) and indices m,n are 
suppressed: 

V • /P = r-2(r2/)'y, V • /B = -r-i/Ay, V • /C = 0; 

A (/A) = (A/) A + r~V^^A where A = P, B, C and (9) 
D^p = 2AB - (A^ + 2) P, D^B = 2AP - A^B, D^C = -A^C. 



4 The Stroh formalism in spherical coordinates 

We consider time harmonic motion with the time dependence e~*'^* omitted but understood. It 
is useful to re-write the equilibrium equation ([1]) in the form distinguishing terms with radial and 
angular derivatives, namely: 

{r'^tr)^r + rr = -r'^puj'^u, (10a) 
with T = (sin 6)^'^ [(sin 9te)^g + t^,<^ + sin 9Ktg + Ht^] , (10b) 

where tj are defined by (j2ap . In the following, the density and elastic coefficients are assumed to 
be radially inhomogeneous. Suppress the indices m, n, denote A (e^) = P, B, C and let 

u = J2j^UaA, tr = Y,j^TAA, T = J^^rAA; (11a) 
V = {Up,UB,Ucf , T = {Tp,T^,Tcf , T = (Fp, Tb, Tc)^ . (lib) 

Our objective in this paper is two-fold: first to find the most general symmetry which admits of 
the separation of variables solution u = Ua {r) A; and second, to obtain these separation of 
variables solutions when they are possible. 

In view of eqs. (|10p - (|lip . the first problem implies answering the following question: given 
the ansatz U = U (r), what symmetry yields simultaneous conditions T = T (r) and F = F (r)? 
Direct calculation of Ta = • A from ()2ap with u = Ua (r) A shows that T = T (r) holds for 
tetragonal, cubic and transversely- isotropic symmetries if their principal axes are parallel to e,., but 
it is invalid for any other cases including the above symmetries with non-radial principal axes, the 
trigonal symmetry and certainly any lower symmetries. Note that the other traction vectors tg and 

do not admit a similar expansion in A with scalar coefficients depending on r even if the material 
is isotropic and uniform; however, this is not so for their combination r. Calculation of Ta = t ■ A 
shows that F = F (r) is possible, but it holds only for transverse isotropy with the principal axis 
along e^.. Thus, altogether the separation of variables solution u = "^j^Ua (^) A is ensured only 
in the presence of transverse isotropy about e^; otherwise the coefficients Ua must depend on the 
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spherical angles (apart from the theoretical possibility of coincidental equalities between material 
constants that is of no practical interest). 

Regarding the second problem, it is solved by deriving the first-order ordinary differential system 
in r (solutions of such ODS are standard). This is done in detail in the next section. For now we 
simply formulate the overall result. 

Theorem 1 Time harmonic solutions of the equations of linear elasticity in a spherically anisotropic 
radially inhomogeneous body with p (r) and Cijki (r) admit of separation of variables using the 
vector spherical harmonic functions only if the material is transversely isotropic about e^. In 
this case, the separable solution for the displacement and radial-traction vectors is in the form 
u = X]a^a('')A, tr = X]a-^a(?')A (with A = P, B, C and the indices m, n omitted). The 
amplitudes U (r) = {Up,Ub,Uc)'^ and T (r) = {Tp,T-B,Tc)^ are defined by the Stroh-like ODS 

V' = -^n with n (r) = ^^^2T J ' ^ = (Q - nT~"R^) - r'pu:^! -^rnT~' ) ' (''"^ 
where ' = d/dr, 

/ 2ci2 Xcqq \ 

T = diag[cii, C66, cee], 7^ = -Aci2 -cee , (12b) 

V -C66/ 

(A^C66 + 4 (C22 - C44) A(2C44 - 2C22 - Cm) \ 

A(2C44 - 2C22 - cm) >?C22 + cm - 2C44 =2^' (^^^^ 

A2c44 + C66 - 2C44/ 

and A = [n{n + 1)]^^^. For real material parameters (and real uP'), G = TG"'"T where ^ means 
Hermitian conjugation and T is a matrix with zero diagonal and identity off-diagonal 3x3 blocks. 



5 Derivation of the Stroh-like ODS for spherical transverse isotropy 
5.1 Elastic coefRcients and stress 

Having established that the separation of variables solution via the use of harmonics works only 
for transverse isotropy, our purpose here is to obtain the explicit form of the ODS which is stated 
by Theorem [H The derivation proceeds by splitting the problem into the isotropic and anisotropic 
parts based on the standard form of the transversely isotropic elastic coefficients Cijki (r) |19j : 

(iso) , (anis) (iso) j- j- , /re- , jf r \ 

4jkT^ = "1 i^ik^jl^ll + ^il^jlhl + ^jk^il^ll + SjAlhl) (13) 

+ 0,2 Sii6ji6kiSii + 03 {6ij6kiSii + dkidii6ji) with 
ai = Cm - C44, a2 = cn + C22 - 2ci2 - 4c66, as = C12 - C23, C22 = C23 + 2C44. 
This partitions stress as 

a = a'^^'°^ + a^^'^''^) , = C23 (V • u) I + 2c44e, 

^(anis) ^ (er 7^ + 7.^ er) + a2Ur,rer (8) e^ + 03 [tir-,rl + (V • u) e^ Cg" er] (14) 

where 7^ = e^e = - (u^^ + Vur — r~^u±) with u_|_ = ugeg + u^e^. 

Another useful split occurs due to uncoupling of the shear wave motion with u polarized along C. 
We begin with this observation and subsequently examine solutions with u spanned by P and B. 
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5.2 Uncoupled SH solution 

Keeping the conventional term 'shear- horizontal' (SH) for the waves with u parallel to C, assume 
that 

nsH = Uc{r) C. (15) 



Inserting the ansatz (|15p into ()14p readily determines the SH radial traction as 

t^ = TcCwith Tc = cqq{U'c - r-^Uc) (16) 

(note that the SH tractions and are not aligned with C). Applying the divergence identities 
([9]) to the isotropic and anisotropic parts of the stress gives the equalities 

diwaf°^ = [<^i^{U'c - r~'Uc) + cm{AUc - X^r~^Uc)] C, (17) 
diYCT^g^P = [a[{Ui. - r-^Uc)C + ai{AUc - 2r-^Uc)] C. 

Adding them leads to the equation of SH motion div ctsh = —p^'^^SH in the form 

{r^ceeU^y + {r^pu^ - X^a^ + 2(c44 - cge) - rc^6)f/c = 0. (18) 

The latter can be recast, using (fT6]l . as 

[r^c^U'c]' - (2C66 + r46)?7c = (r^Tc)' + rT,. (19) 

Thus (jl6p is the first and ()19p the second equations of the following ODS for SH waves: 



Uc\ _■( -ir~^ -(r%6)-^\ ( Uc 

ir^TcJ \{X^ - 2)c44 - r^poJ^ ir'^ ) [ir^Tc 



This is clearly the same as the pair of equations in the third and sixth rows of system ()12ap . 
5.3 In-plane problem 

We now consider U = UaA. and T = 7a A where A = P, B only, in which sense this case 
may be referred to as the 'in-plane' problem. It implies restricting attention to the upper 2x2 block 
of the matrices in (11211 . 



Isotropic part Assume 

^in-plane = J^A ^^^^^^ A = P, B. (21) 

The equation for t^ = erff readily follows from o"('^°^ in (jl4p i to yield the isotropic part of the first 
equation of the in-plane ODS as 

.T,„ = rT^V' + wHh ^„ = ( - ,° J . = ( f - - ) ■ 

In order to formulate the second equation, which is based on the equilibrium condition div<T = 
— pw^u, it is convenient to use the following identity for the radially inhomogeneous medium: 



t,) +rVu= pt,) -r^diva , (23) 

' L ' J homo 
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where the 'homogeneous' term on the right hand side is understood to be evaluated as if the 
elastic moduli are independent of r, i.e., derivatives of Cijki{r) are ignored. This reduces the task 
to concentrating on the right member in (j23p . namely, to expanding divcThomo in vector spherical 
harmonics given that u = J^aA: 

divffhomo = J^^i^AA withA = P,B. (24) 

Using divo"[j'^^Q = (C23 + C44) V (V • u) + C44V^u and the identities ([9]) yields the isotropic part Fiso 
of the vector F = (Fp, Fb)^ as 

_ /C22A - ^ (2C22 + A2C44) ^ [(C44 " C22) ^ + \ (C22 + C44)] \ 

" U [(^22 - C44) i + ^C22] C44 - gc22 ) ^^^^ 

Combining ()25p with (122^ and substituting into (j23p gives 

(r^Ti^o)' = r7^i3oU' + (Q.o - rVo^^l) U with Q,.^ = " ^^^^ + ^'^^^ ^ (^44 - 2c22)\ _ 

^ ^ ^ ^ V ^ (C44 - 2C22) -C44 + A^C22 / 

(26) 



Anisotropic part From and the definition of Oj in (I13p . the anisotropic part to be added 

to the first equation is 

rTanis = rTanisU' + 7^LsU with (27) 

^ ^ /^Cii - C23 \ ^ 2 (Ci2 - C23) A (C23 - 012] 

\ C66 - C44y ' \A (C66 - C44) C44 - C66 



Again using (fM|l '?. now with the identities 

2 d 2 2 

V [(V • u) Qr^Gr] = (A + j)?!, ^e^, V [/(r)e^®er] = (/,^ + -f)er, (28) 

V [/ (r) (e^0B + B ® e^)] = /'B + /V (e^«)B + B (g) e^) = /'B + r"!/ (3B-AP) , 
yields the anisotropic part Fanis of the vector F = (Fp, Fb)"^ defined in (j24|) : 

_ /(a2 + 2a3) A + (2a3 - 4ai - AV) ^ [(ai - 03)^: + i (3ai - 03)] \ . . 

■c anis — I A// \d,2\ I \ 2\ \ \'^^) 

\ 7((a3-ai)df + ai(A-;7) y 



Combining (j29j) with (p7j) and substituting the result into (|23j) gives 



(r^Tanis)' = r7^anisU' + (Qanis " r^^^l) U with Qanis = (c66 " C44) (^^^ • (30) 



Result: The isotropic and anisotropic parts of the in-plane solution can now be superimposed. 
Adding ([27]) to ([221) and ([30]) to (l26|) gives the in-plane part of (fT2all . This, together with ([20]) . 
completes the proof of Theorem [TJ 
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6 Discussion 



6.1 Radially uniform materials 

Consider first the SH motion. Assuming radially uniform material coefficients and denoting Uc{r) — s- 

1 /2 

u{R) with R = ruj^p/cee) , the uncoupled equation p8|) becomes 

(R^u'V + [R^ - fi{fi + l)]u = 0, + = [n(n + l) -2]— + 2, (31) 

C66 

the solutions of which are spherical Bessel functions j^{R) and y^iiR) (alternatively, Hankel func- 
tions /ii^^(i?), h!i?{R)). The identity ([311)2 was obtained by [8] for static equilibrium, and the vector 
function defined by (llSp with Bessel or Hankel function solutions for Uc{r) is equivalent to the 
vector function M{r,6,(j)) of [18] (see also [20]). 

The in-plane differential equations for displacement and traction amplitudes may be recast as 
second-order differential equations for Up and Ub- In the case of radially uniform materials, these 
equations reduce to 



/ puj'^r'^ - 2(c22 + C23 - C12) - A^cee A(c22 + C23 - C12 + cm) \ fUp\ (^2) 

V A(C22 + C23 + 2C66) pUJ^r^ + 2(c44 " Cqg) - X^C22 J \Ub J ' ^ ' 

Even though the material coefficients are constant, this coupled pair of equations does not seem 
to admit an evident explicit solution in terms of special functions and should be solved by means 
discussed in ^6.2p . For isotropic materials with cj^ = cn/p and c|. = cqq/p, eqs. ([32]) further 
simplify to the form 

^2 

{r^WpY + i—r^ - - 2)Up + 2XUb - Xcl^{cl - c^tWUb)' - AC/p] = 0, (33a) 



^2 1 
{r^U's)' + (— - X^)Ub + 2XUp + Xc^^cl - 4)[-(r2[/p)' - XUp] = 0, (33b) 

Crp f 



which leads to the known isotropic solutions. They may be elicited as follows. Set the final term 
in (|33ap to zero by assuming for some as yet unknown function v that 

Up = v\ Ub = X-. (34) 
r 

Substituting this into equations ([33]) recasts them, respectively, as 

d 2 A ^2 

(- )Mlv = 0, — ^Mlv = 0, with Mlv = (r^v'y + (^r^ - X'^)v. (35) 

M r r ^ r cj. cj^ 

Both equalities are satisfied if i; is a solution of the same spherical Bessel equation, M^v = 0, which 
has solutions jn{ruj/cL) and yniruj/cL,)- Similarly, setting to zero the final term in ([33bp by taking 

w 1 

Up = A-, Ub = -irw)', (36) 

for some function w, we find that ([33aD and (p3bl) b ecome, respectively, 

\„2 d 1 up' 

— "^Mtw = 0, (- )Mtw = 0, with Mtw = (r'^w')' + (^r^ - X^)w. (37) 
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These are satisfied if ti; is a solution of the spherical Bessel equation of order n, Mtw = 0. 

The uncoupled longitudinal and transverse wave solutions ()34p and ()36p for the uniform and 
isotropic case are consistent with the potential representation using the Helmholtz decomposition. 
The vector functions of r,6,(j) formed from (f2T|) with (p4|l . (f36ll can be identified as the vector 
functions L and N, respectively, of [18] (see also |20j). 



6.2 Radially inhomogeneous materials 

The first-order ODS (jl2ap with material coefficients depending on r admits a general solution in 
the standard form r/(r) = M(r, ro)r/(ro), where 7?(ro) is the initial data and the matricant M(r, rg) 
may be evaluated by the Peano series [21] 



Mir, ro) = 1+ / d.^ + r dx^ f dx,^ + ... (38) 



The matricant solution applies for r, ro 7^ 0. The case where the solution needs to be extended to 
the origin point r = requires a special treatment based on the theory of ODS with an irregular 
singular point [22]. A similar state of affairs arises for the Stroh-like ODS in cylindrical coordinates 
at the axial points r = (except that they are regular singular points) [23^ 17]. In the following we 
assume r, ro / 0. 

By analogy with the cylindrical case, the algebraic symmetry G = TG~''T (see Theorem [1]) 
yields M+TM = T and 

^ [M+TM) = with M{r) = , (39) 

where the 6x6 matrix M{r) is composed of any six linearly independent solutions ?7^"-*(r) = 
(UW, ir^T^")) , a = 1..6, of ([US]). Thus AA+(r)T7V(r) is a constant matrix which can be chosen 
to provide the partial solutions ?7^"^ (r) with appropriate pointwise orthogonality in the sense of a 
product ([MD [6]- 

Recalling the angular indices of the vector spherical harmonics Amn(er.)=Pmni Bmn> Cmn 
defines the displacement-traction modes in full as 



(a) \ °° 

^^"^«=( ]=Y.^^n^ (a = 1...6) with 

I n=0 



where it is taken into account that U^'^ ^, T^'^ ^(r) are independent of m (see ^6tl6.3p ). There is 
no pointwise orthogonality of the 'full' modes in general due to the complex conjugation in (j39p and 
its absence in ([7]). At the same time, one can make use of ()39p along with the integral orthogonality 
of vector harmonics ([8]). This enables evaluation of the angular (and time-period) average of the 
radial component of energy flux associated with the mode rf^"'^ (r) as 

P^°'^ = - j ^t(°) •u(")^^d(7 = j V°)+(r)Tj7(")(r)dO. (41) 

This may be further reduced by using Q, (139]) and (f4T]) . 
00 „ 
^'J"^ = E^- ' whereP(^) = - / (t^ • u^f)) dl7 = (2n + l),7^r^+(r)Tr7(f)(r). (42) 

n=0 ■'^ 
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Thus the fluxes P^r carried by the modes ryn (r) add up to give Pr of glD. According to ([39]) . 
these fluxes do not depend on r which is consistent with energy conservation for the assumed case 
of real material parameters. 



6.3 The case u(r, 6) (m = 0) 

The separation of variables using vector spherical harmonics involves a single auxiliary equation ([5]) 
which does not depend on m, and hence Theorem [1] holds for any m [9]. This does not, however, 
preclude the possibility that some lower symmetry permits such separation for the specific case of 
m = 0. The derivation along the lines of ^ shows that 'radially tetragonal' symmetry, with e,. 
parallel to the 4-fold axis and with C24, C34 = 0, admits separation of variables for m = 0. The 
result for this tetragonal symmetry at m = amounts to replacing 2C44 by C22 — C23 in all entries of 
Q in (jl2cp . except the term A^C44 in the last diagonal component. Any other tetragonal symmetry 
(including that with || 4- fold axis but with C24 = — C34 7^ 0), any trigonal symmetry and indeed 
all lower symmetries prevent separation of variables even for m = 0. 

The same question may be raised concerning the SH uncoupling for m = 0. The above consid- 
eration implies that (j20p with m = can be extended to 'radially tetragonal' symmetry. Moreover, 
it can be shown that orthorhombic symmetry aligned with the coordinate planes also leads to SH 
uncoupling for m = 0, with the equation of motion obtained from (j20p by replacing cge with C55. 
In this case the SH modes admit the separation of variables form \isH{f,0) = Uc{r)X~^ Pn,e{(^)^^ 
while the modes u(r, 6) polarized in the plane {e^, eg} do not. SH uncoupling is generally precluded 
for tetragonal symmetry with C24 = — C34 7^ and for trigonal symmetry even if m = 0. 



6.4 Solutions for n = 

Consider the term t/q"^ (r) with n = (and hence m = 0) of the series ()40p . The corresponding 
spherical harmonics are 

Poo = e„ Boo = 0, Coo = 0. (43) 

Vanishing of Boo and Coo rnay be formally inferred from the definition ^ by taking the limit of 

as — )• 0. It is also consistent with another framework that defines Bmn 
and Crnn without a normalization factor A~^. Note however that the symmetry G = TG+T of 
the system matrix in (I12p is not preserved under any change of the definition ([6]) unless it implies 
multiplying all three harmonics A^n by the same factor. 

By ()43p . the sextet t/q"^ (a = 1..6) for n = contains only two non-zero modes which are the 
modes a = 11,12 of the longitudinal (radially polarized) wave 

u = U{r) e,., = r(r) e^, (44) 

where we drop the indices used in (|40p . Equation (jl2ap with n = provides an uncoupled system 
of two equations for the amplitudes U and T: 



U \' . ( 2ir ^ci2q/ 



ir'^T) ' \. 4 ( C22 - C44 - ^ ) - r'^poj'^ -2ir ^ci2C 




(45) 



which reduces to a single second-order equation for f/(r), 

{r^ciiU')' + {r^poj'^ + 4(c44 - C22) + 2(rci2)')f^ = 0. (46) 
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Note that the wave (j44p exerts the tractions tj = f{r)ei for i = 9,(j), where / = C12U' + 
(c22 — C44) U. Equations (HSl) and ([I6]l are similar to the SH-wave eqs. (f20l) and (fT9]) formahy 
taken with n = 0. Since the wave ()44p also implies m = 0, it follows from §6(c) that eqs. (j45p . (|46p 
obtained from the transversely isotropic system (|12bp and (jl2cp can be extended to the tetragonal 
case, by replacing 2C44 with C22 — C23. The same results can be obtained by inserting u = U{r)er 
in the initial elastodynamic equations ([I])-© with tetragonal symmetry. Note in this regard that 
the elastostatic solution of the form ()44p was analyzed for various symmetries in \13\ I24| and |25| . 

Interestingly, replacing in the harmonics definition ([6]) with some other power of A may 
lead to peculiar solutions at n = different from (|44p . This is a problem of interest in its own 
right. Consider the set of spherical harmonics defined as 

The corresponding spherical harmonic functions for n = are (dropping the m = subscript) 

0-0 1 
Po = AoBr, Bo = -^0 tan - e^, Co = ^0 tan - e^, with ^o = 7—;=, (48) 

which follow from (j47p using first an exchange of limits, 

1 d „ , d , P^(cos6i) - Po(cos6') ,^ , 

lim — TT^Pi. cos 0) = — lim ^ / ^ ^, 49 

and then Jolliffe's formula [26] to evaluate the derivative with respect to 



F„(.) = -P„(.) 1„ i±i + ^i^, - D" 1„ i±i ) . (50) 



By (I48p . Po corresponds to purely radial motion, while Bo and Co (non-zero in contrast to (I43p ) 
represent shearing and twisting about the polar axis, with zero at one pole and a singularity at 
the other. The latter can be avoided by introducing a conical cut centered at r = of arbitrarily 
small angular extent at ^ = tt. This obviates the singularity allowing the otherwise un-normalized 
spherical harmonics Bq and Cq to satisfy the orthonormality conditions ([8]) . Note that the spherical 
harmonics Bq and Cq are derived from and related with Legendre polynomials of the first kind 
which define the spherical harmonic Pq, although the Legendre polynomials of the second kind may 
be represented as Qn{z) = \\Fn{.z) - (-l)"F„(-z)] (see eq. ([MD) [SH]. 

Based on the above, it is of interest to consider specific n = solutions in the form 

^ = Eao^Ao(OAo, t. = ^^^7X^,(r)Ao, Ao = Po,Bo,Co. (51) 

The system (jl2p modified with respect to the spherical harmonics A of (j47p and taken for n = 
defines the displacement amplitudes C/^^ (r) by the following equations: 

(r^cnC/^J + [r'pu? + 4(c44 - C22) + 2{TCx2)')Up, = [2(044 - C22 + — ) - rc^jg] U^^ 

^11 

+ r2[r-i(c66 + ci2)?7Bj', (52a) 
(r'c66?7^j' + (rV^' + 2(C44 - C66) " rc'ge) U^^^ =0, Aq = Bq, Cq. (52b) 

It is seen that eqs. (j52bp for the amplitudes U-q^ and are the same as the SH wave equation 
(fTSl) with n = 0, and that C/g^ provides a forcing term in eq. (I52ap for Up^ (note also that (152aP 
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reduces to ([16|) if U^^^ = 0). Thus the SH solution u = U^^^Cq (|| e^) is completely uncoupled. 
However, the SH motion polarized in the Bq || eg direction drives the radial motion and this results 
in the coupled wave u = Up^Po + J/g^^Bo- At the same time, non-zero angular motion requires that 
the issue of singularities at ^ = vr be resolved. Under normal circumstances, e.g. a solid sphere, or 
a complete shell, this is not the case, and any angular motion is precluded, i.e. U-^^^ = = 
and only radial motion u = U (r) e,. occurs (where U is defined by (|46p ). If the small conical 
cut device is introduced, then the surface of the cone must apply a force and a moment sufficient 
to maintain the dynamic tractions required of the solutions. The magnitudes of the tractions for 
different types of motion depend critically on the cone angle e <C 1. For purely longitudinal motion 
the traction is independent of e, i.e. 0(1). The (e,., e0)-coupled motion requires normal traction 
tg of 0( tan^ f = 0(e~^) and is therefore ruled out as a viable n = dynamic solution. 

In the case of the pure SH motion u = U^^Cq, the normal traction tg = r~^C44 tan | u is 0(e~^), 
a pure twist in the direction corresponding to a net torque about the polar axis of order unity 
(re X 27rre x This suggests that pure twisting motion may be induced in a solid sphere with 

a fixed polar axis by application of torque to a small conical insert. 

7 Conclusion 

The central result of the paper, Theorem [H shows that spherically anisotropic radially inhomo- 
geneous materials admit elastodynamic solutions u(r, 9, <j)) in a separation of variables form. The 
angular dependence (on 0, (j)) is described by the vector spherical harmonics while the radial depen- 
dence (on r) is separated and determined by the Stroh-like first-order ordinary differential system 
that is solvable by standard means. It is proved that such separation of variables solution is gener- 
ally possible only if the spherical anisotropy is restricted to transverse isotropy with the principal 
axis in the radial direction e^. Transverse isotropy about distinguishes a class of materials with 
complete spherical anisotropy, which is a physically natural model of spherical anisotropy since it 
ensures invariance of material properties with respect to any orientation of the polar axis. 

The separable of variables solution u(r, 0, 4>) for the transverse isotropic case does not explicitly 
depend on the azimuthal order m. At the same time, dependence on m reveals itself in that the 
solution of the form u(r, 6) (i.e. with m = 0) admits separation of variables via spherical harmonics 
not only for transverse isotropy but also for lower symmetry - but only 'up to' tetragonal with 
along the 4- fold axis and with C24, C34 = 0. Note that the solutions u(r, 6, (p) for transverse isotropy 
and u(r, 9) for the above tetragonal symmetry uncouple the shear modes parallel to the vector 
harmonic C (|| e,^ at m = 0) from the in-plane modes orthogonal to C. Moreover, shear modes of 
the form u(r, 9) are uncoupled also for orthorhombic symmetry but not for trigonal symmetry nor 
for tetragonal if C24, C34 / 0. 

The establishment of the Stroh format for the elastodynamic equations in spherical coordinates 
opens the door for applications to various boundary value and scattering problems. For instance, 
solutions for acoustic and elastic wave scattering from solid spheres and shells, which have been 
limited to isotropic materials (see \17\ §4.10] for a review) or transversely isotropic shells with 
m = [TTj, can be generated for arbitrarily layered shells and solids using standard solution 
techniques outlined in ^6l !j6.2p . Other possible approaches that can be explored with the Stroh 
formalism include impedance matrices for spherical shells and solids, the use of which simplifies 
the formulation of boundary value problems, e.g. determining modal frequencies, solving radiation 
and scattering problems. By analogy with the cylindrical situation [7j it should be possible to 
formulate a matrix Riccati ordinary differential equation for the impedance matrix as a function of 
the spherical radius r, with a unique solution at the origin that depends only on the elastic constants 
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at r = 0. The Stroh formalism is particularly suited to solution of elastodynamic problems with 
forcing, e.g. from thermal expansion via laser excitation with application to non-destructive testing. 
More exotic issues could be addressed with the Stroh system, such as modelling and simulation 
of fully elastic 3D 'radial wave crystals' [27], i.e. shells of radially periodic materials that exhibit 
Bloch wave effects normally associated with rectangular periodic crystals. 

Acknowledgements 

A.N.N, acknowledges support from CNRS and from ONR. 



References 

[1] J. J. Wu, T. C. T. Ting, and D. M. Barnett, editors. Modern Theory of Anisotropic Elasticity 
and Applications. SIAM, Philadelphia, 1991. 

[2] T. C. T. Ting. Anisotropic Elasticity: Theory and Applications. Oxford University Press, 
1996. 

[3] K.-C. Wu. Generalization of the Stroh formalism to 3-dimensional anisotropic elasticity. J. 
Elasticity, 51:213-225, 1998. 

[4] A. L. Shuvalov. On the theory of wave propagation in anisotropic plates. Proc. R. Soc. A, 
456:2197-2222, 2000. 

[5] A. L. Shuvalov, O. Poncelet, and M. Deschamps. General formalism for plane guided waves 
in transversely inhomogeneous anisotropic plates. Wave Motion, 40(4):413-426, 2004. 

[6] A. L. Shuvalov. A sextic formalism for three-dimensional elastodynamics of cylindrically 
anisotropic radially inhomogeneous materials. Proc. R. Soc. A, 459(2035):1611-1639, 2003. 

[7] A. N. Norris and A. L. Shuvalov. Wave impedance matrices for cylindrically anisotropic radially 
inhomogeneous elastic materials. Q. J. Meek. Appl. Math., 63:1-35, 2010. 

[8] H. C. Hu. On the general theory of elasticity for a spherically isotropic medium. Acta Sci. 
Sin., 3:247-260, 1954. 

[9] N. A. Shul'ga, A. Y. Grigorenko, and T. L. Efimova. Free non-axisymmetric oscillations of 
a thick- walled, nonhomogeneous, transversally isotropic, hollow sphere. Int. Appl. Mech., 
24(5):439-444, 1988. 

[10] W. Chen and H. J. Ding. Free vibration of multi-layered spherically isotropic hollow spheres. 
Int. J. Mech. Sci., 43(3):667-680, 2001. 

[11] S. Hasheminejad and M. Maleki. Acoustic wave interaction with a laminated transversely 
isotropic spherical shell with imperfect bonding. Arch. Appl. Mech. (Ingenieur Archiv), 
79(2):97-112, 2009. 

[12] C. Scandrett. Scattering and active acoustic control from a submerged spherical shell. J. 
Acoust. Soc. Am., 111:893-907, 2002. 

[13] T. C. T. Ting. The remarkable nature of radially symmetric deformation of spherically uniform 
linear anisotropic elastic solids. J. Elasticity, 53:47-64, 1998. 



14 



[14] S. G. Lekhnitskii. Theory of Elasticity of an Anisotropic Elastic Body. Holden-Day, San 
Francisco, 1963. 

[15] S. Nair. Introduction to Continuum Mechanics. Cambridge University Press, New York, 2009. 

[16] J. Lothe and D. M. Barnett. On the existence of surface-wave solutions for anisotropic elastic 
half-spaces with free surface. J. Appl. Phys., 47(2):428-433, 1976. 

[17] P. A. Martin. Multiple Scattering: Interaction of Time-harmonic Waves with N Obstacles. 
Cambridge University Press, New York, 2006. 

[18] P. M. Morse and H. Feshbach. Methods of Theoretical Physics, Vol. II. McGraw-Hill, New 
York, 1953. 

[19] F. I. Fedorov. Theory of Elastic Waves in Crystals. Plenum Press, New York, 1968. 

[20] G. Dassios and Z. Rigou. Elastic Herglotz functions. SIAM J. Appl. Math., 55(5):1345-1361, 
1995. 

[21] M. C. Pease. Methods of Matrix Algebra. Academic Press, New York, 1965. 

[22] W. Wasow. Asymptotic Expansions for Ordinary Differential Equations. Interscience, New 
York, 1965. 

[23] A. L. Shuvalov. The Frobenius power series solution for cylindrically anisotropic radially 
inhomogeneous elastic materials. Q. J. Mech. Appl. Math., 56(3):327-345, 2003. 

[24] S. S. Antman and P. V. Negron-Marrero. The remarkable nature of radially symmetric equi- 
librium states of aeolotropic nonlinearly elastic bodies. J. Elasticity, 18:1317-164, 1987. 

[25] S. S. Antman and T. C. T. Ting. Anisotropy consistent with spherical symmetry in continuum 
mechanics. J. Elasticity, 62(l):85-93, 2001. 

[26] A. E. Jolliffe. A form for ■^Pn{fJ,), where PnifJ') is the Legendre polynomial of degree n. Mess. 
Math., 49:125-127, 1919. 

[27] D. Torrent and J. Sanchez-Dehesa. Radial wave crystals: Radially periodic structures from 
anisotropic metamaterials for engineering acoustic or electromagnetic waves. Phys. Rev. Lett., 
103(6):064301-h, 2009. 



15 



